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Abstract: A quantum harmonic oscillator (spring subsystem) is stabilized towards a target Fock 
state by reservoir engineering. This passive and open-loop stabilization works by consecutive and 
identical Hamiltonian interactions with auxiliary systems, here three-level atoms (the auxiliary 
ladder subsystem), followed by a partial trace over these auxiliary atoms. A scalar control input 
governs the interaction, defining which atomic transition in the ladder subsystem is in resonance 
with the spring subsystem. We use it to build a time-varying interaction with individual 
atoms, that combines three non-commuting steps. We show that the resulting reservoir robustly 
stabilizes any initial spring state distributed between and 4n + 3 quanta of vibrations towards 
a pure target Fock state of vibration number n. The convergence proof relies on the construction 
of a strict Lyapunov function for the Kraus map induced by this reservoir setting on the spring 
subsystem. Simulations with realistic parameters corresponding to the quantum electrodynamics 
setup at Ecole Normale Superieure further illustrate the robustness of the method. 
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1. INTRODUCTION 

The last decades have seen a surge of developments on 
the control of systems that feature essential quantum dy- 
namics (see e.g. Wiseman and Milburn [2009]). A major 
motivation is to harness their peculiar possibilities for IT 
applications ranging from fundamentally-secure communi- 
cation, over hyper-precise measurements, to the quantum 
computer (see e.g. Nielsen and Chuang [2000], Haroche 
and Raimond [2006] ) . A basic building block for operations 
with quantum dynamics is the ability to produce and 
stabilize a wealth of different 'target' states. The extreme 
fragility, limited measurement and control possibilities in 
quantum systems make this already a challenging task for 
many target states with non-trivial (and thus interesting) 
properties. 

Since an isolated quantum system fundamentally follows 
pure Hamiltonian dynamics, asymptotic stabilization nec- 
essarily requires interaction with the external world and 
relies on the theory of open quantum systems. Such in- 
teractions can involve measurement and feedback action, 
called measurement-based feedback, or also control by tai- 
lored interaction as broadly advocated in Willems [1995] 
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and called coherent feedback in the quantum context. 
Measurement-based feedback introduces specific difficul- 
ties because quantum measurement is fundamentally lim- 
ited to partial state knowledge and always perturbs the 
measured system ( "back-action" ) . Controller design there- 
fore always needs to follow an interactions-based reasoning 
(see e.g. Dotsenko et al. [2009], Sayrin and et al. [2011]). 
Coherent feedback involves a specific structure, based on 
joint evolution on a tensor product of the interacting sub- 
systems that leads to dissipative and/or stochastic quan- 
tum dynamics for the target subsystem (see e.g. James 
and Gough [2010], Kerckhoff et al. [2010]. 

Reservoir engineering is a systematic method related to 
coherent feedback for open loop stabilization of quantum 
systems through tailored interaction. For discrete-time 
quantum systems, the target system repeats the same 
interaction with a succession of auxiliary systems, which 
are discarded after interaction. The dynamics on the tar- 
get system resulting from repeating the same interaction 
process, takes the form of a Kraus map, see Kraus [1983]. 
A random Kraus map would generally drive the system 
to a mixed quantum state (featuring classical uncertain- 
ties) . For reservoir action, interaction with each individual 
auxiliary system must be tailored such that the Kraus 
map has a desired pure state as asymptotically stable 
equilibrium. General techniques have been proposed to 
tailor Kraus maps, provided any chosen unitary evolutions 
can be applied to the target subsystem by controlling its 
Hamiltonian (kind of 'complete actuation' on the sub- 



system Hamiltonian, see e.g. Ticozzi and Viola [2009], 
Bolognani and Ticozzi [2010]). There are however many 
systems where control over the Hamiltonian is instead very 
restricted. 

This paper considers such a situation, where the tar- 
get system is an electromagnetic field mode, interacting 
with three-level atoms as auxiliary systems. This quan- 
tum electrodynamics situation is a prototype for so-called 
universal spring-spin systems (see e.g. Haroche and Rai- 
mond [2006]). The spin-spring Hamiltonian is controlled 
by shifting the frequencies of the atomic transitions as 
a function of time, through the Stark effect. This yields 
a very low-dimensional input signal to tailor the target 
system's evolution governed by the resulting Kraus map. 

Sarlette et al. [2011] propose a method to stabilize so- 
called 'Schrodinger cat states' with this setup. It builds 
the overall unitary interaction operator as a symmetric 
product of non-commuting basic operators, by varying the 
input signal during the field's interaction with each atom. 
We prove here that stabilizing pure photon states ('Fock 
states') is also possible. Rempe et al. [1990] have proposed 
a method based on 'trapping state conditions', where the 
auxiliary systems are two-level atoms interacting reso- 
nantly with the quantized field mode: the target Fock state 
is an equilibrium of the Kraus map but it is a saddle 
point with stable manifold in some but not all relevant 
directions. The present paper proposes a modification of 
this method that makes the target Fock state a stable 
equilibrium of the Kraus map in the most relevant Hilbcrt 
subspace. We therefore build an interaction that combines 
the 'trapping state' approach of Rempe et al. [1990] with 
the construction of symmetric products of non-commuting 
operators from Sarlette et al. [2011]. 

We give formal convergence properties for the resulting 
scheme, both in absence and in presence of a disturbing 
environment (Theorem 1 based on the construction of a 
strict Lyapunov function, and Proposition 2). In section 2 
we define the target system, the auxiliary three-level 
atoms, the interaction Hamiltonian with its scalar control 
input u and the associated Kraus map. In section 3, the 
operators appearing in the Kraus map are computed for 
our open-loop piecewise constant control (10). Section 4 is 
devoted to convergence analysis. Simulations in Section 5 
briefly explore the influence of parameter uncertainties and 
illustrate the robustness for realistic QED parameters. 

2. SYSTEM DESCRIPTION 

We consider a quantum electrodynamics setup as de- 
scribed e.g. in Haroche and Raimond [2006]. The target 
system is a field mode, with infinite-dimensional Hilbcrt 
space H = span(|0), . . . , \n), ... ) and free Hamiltonian 

H f = ^2 n Uf\n)(n\ (1) 

71 = 

where Uf > is the field mode pulsation and the 
index n gives the photon number. This paper denotes 
by \n) the pure photon number cigenstates, called Fock 
states. The photon number operator is defined by N = 
J2n=o 71 \ n ) ( n l = Hf I w f- W e will further use the notation 
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/n = /(N) = ]T f(n) \n)(n\ = £ /„ \n)(n\ 

n=0 n=0 

for an operator that is diagonal in the Fock basis, for any 
function 1 / : N ^ C. Wc denote by U™\ the Hilbert 
subspace spanned by Fock states |ni),|ni + 1), \n-i). 
A stream of identical atoms consecutively interact with 
this field mode according to the Jaynes-Cummings model, 
playing the role of a 'reservoir' to stabilize the field in a 
target state. We consider three atomic levels \g), |e), |m), 
with free Hamiltonian 

H a = uj m \m)(m\ - ui g \g){g\ (2) 

up to an irrelevant term that is a multiple of the iden- 
tity operator. The transition frequencies between levels 
(\g), |e)) and (|e), |m)), that is uj g <G R and oj m € R respec- 
tively, depend on an input u € R. The latter represents an 
energy shift induced on |e) by an external field through a 
Stark effect, such that 

0Jg(u) =U g + u and oj m (u) = u] m — u . (3) 

The constants uJ g and uJ m are the transition pulsations 
in absence of external field. For simplicity, the following 
assumes that Lu g (u) > and u> m (u) > 0: \g) corresponds 
to the lowest atomic level whereas \e) and |ra) are the 
first and second excited states (3-level ladder system). 
The proposed method can be adapted to other energy 
arrangements, like V-structures, by adding short external 
control pulses that switch the atomic state at specific 
points in our scheme. 

To make atomic transitions and field mode interact, we 
consider the following realistic situation, in the spirit of 
Santos and Carvalho [2011]: \ZJ g — \uJ m — w/| -C w/; 
\u\ <~ \UJ g — 0Jf\, \u\ ~ \uJ m — uif\. Then with the standard 
rotating wave approximation, the atom-field interaction is 
described by the Hamiltonian: 

H c = i ^ (at (\g)(e\ + \e)(m\) - a(|e)( 9 | + |m)(e|) ) . (4) 

Here i = >/— T and f2 is the interaction strength 
factor, assumed to be equal for transitions (g, e) and 
(e, m ); photon annihilation operator a is defined by a = 

J2n=i V™ l n "l)( n l m the Fock basis; and t denotes the 
adjoint of an operator (complex conjugate transpose of the 
associated matrix). We have the fundamental identities 
a^a = N, as well as a /(N) = f(N + I) a and its adjoint 
/(N) at = at /(N + I) 'for any function / : N H C. Equa- 
tion (4) essentially expresses that atomic state can raise 
from \g) to |e) or from |e) to \m) by absorbing one photon 
from the field, or fall inversely by releasing one photon. 
The propagator U, expressing the transformation that the 
joint atom-field state undergoes during interaction, follows 
the Schrodinger equation 

£ t V(t) = -i(H f + H a + H c )V(t) (5) 

with initial condition U(io) = I the identity operator. 
A standard change of variables \if>) — > e lHft \tjj) on field 
state and (\g), |e), \m)) (e~^|g), |e), e^|m)) on 

atomic states leads to 'interaction coordinates'. In these 
coordinates the propagator follows 

&\J(t) = -iHj C U(t) (6) 



1 We take the convention, especially useful in quantum mechanics, 
that N includes 0. 



where the Jaynes-Cummings Hamiltonian writes 
Hj C {u) = (A m - u) \m)(m\ - (A g + u) \g){g\ + H c 



(7) 



with A m = ZJ m — u>f and A g = W g — uif. The goal 
being to allow separate interactions of the two atomic 
transitions with the field, we assume that the |m) and 
\g) levels are attributed sufficiently different frequencies 
in Hjc, i.e. |A fl + A m | =: A » Cl. Note that H JC acts 
in parallel on a set of decoupled subspaces spanned by 
(\g) <g> \n + 1), |e) (g |n), \m) (g |n — 1)). The associated 
matrix operators are thus block-diagonal with blocks of 
size 3x3 at most, which facilitates analysis. 

The atoms are sent one after the other, every T s seconds, 
to undergo the same interaction with the field. We denote 
the field mode density operator just before interacting with 
the (k + l)th atom by p k . The initial state |w at ) e C 3 of 
the atoms can be chosen, as well as the Stark detuning 
signal u(t) during interaction time [0,T], with T < T s . 
Denote Ut the solution at time T of (6) with U(0) = I 
and with the chosen u(t), governing a time- varying Hjc 
in (7). Then the atom-field joint state just after the k + 
1th interaction is given by TJt (pk <g |«at)(u a t|) Ut- This 
generally corresponds to an entangled situation. Since we 
do not measure the final atomic state, the expected field 
evolution follows the Kraus map 

p k+ i = $(p fe ) = Mgp k Ml + M ePk M\ + M mPk Ml , (8) 

where the operators M g , M e , M m acting only on the field 
mode are identified from Ut and |u a t): V|^>) e 'H, 

U T |V)Kt) = M S |V> \g) + M e |V) \e) + M m |V> |m). 
The goal of open-loop stabilization by reservoir engineer- 
ing is to select u(t) and |w a t) such that the dynamics (8) 
asymptotically stabilize a target pure state p k — > Poo = 

|^oo)(V>oo|- 

3. CONTROL DESIGN 

Our objective is to stabilize a given Fock state IV'oo) = \n). 
Rempe et al. [1990] have noted that in absence of the 
\m) level (i.e. |A m — u\ 3> 0), a single atomic transition 
(|g),|e)) in perfect resonant interaction with the field 



(i.e. A s 



0) allows to "trap" states below |n). The 



corresponding propagator (obtained by direct integration 
of (6) which is then 2x2 block-diagonal) writes 



Ur =cos(^) \g)(g\ + cos(^S±I) |e)(e| 



(9) 



-a ^ — \e){g\ 



sin(^F) 



N 



N 



^ f l5)(e| 



where 9 r = t r Q is the interaction strength integrated over 
chosen interaction time t r . Then taking 6 r = 2Tr/y/n + 1, 
the operator in (9) implies no exchange between joint 
state components \g) ® \n + 1) and |e) ® \n). The sub- 
space Hq spanned by Fock states |0), |1), ...\n) then re- 
mains decoupled from the rest of the field Hilbert space 
throughout reservoir action. Now take |u at ) = |e), thus 



M, 



sin( — 5 ) + 

■a T , 



g - M e = cos(^H) and M m = 0. 

Then p k following (8) starting from po with support in 
Hq converges to the trapping state \n)(n\ [Haroche and 
Raimond 2006, page 210]. But if the support of p is 
in 'H|+"i 3 , then p k converges to |4n + 3)(4n + 3| as the 



field gets continuously excited. Therefore any fraction of 
density that is pushed above \n) is lost away to high 
photon numbers. Since perturbations will always induce 
transitions between nearby energy states, this makes the 
method of Rempe et al. [1990] unusable in practice as it 
leaves equilibrium \n) unstable in important directions. 
A robust stabilization method should enlarge the basin of 
attraction of \n)(n\ to include at least any p with support 
in %n\ f° r some n i < n < n 2- 

We achieve such stabilization with m = and n 2 = 4n + 3 
(see Theorem 1) by exploiting the possibility of varying 
u(.) during the interaction. Specifically, we take 



u{t) 




forte [0,(T-t a )/2] 
for t e [(T — t s )/2, (T + t s )/2] (10) 
g forte [(T + t s )/2,T]. 

Fast set-point changes can indeed suitably be experimen- 
tally implemented. The switching time t s and overall in- 
teraction time T will be tuned to optimize operation. We 
still take |w at ) = \e) as initial atomic state. 

The propagator Ut, solution of (5) with u given by 
(10), readily writes Ut = Ui U 2 Ui where Ui = 
exp[-iH JC (-A g ) (T - t 8 )/2] is the solution at t = (T - 
t s )/2 of (5) starting at to = 0, with constant u = — A g ; 
and U2 = exp[— iHjc(A m ) t s ] is the solution at t = t s of 
(5) starting at t = 0, with constant u = A m . We compute 
those operators using quantum Hamiltonian perturbation 
theory on the decoupled subspaces spanned by (\g) (g \n + 
1), |e) ® \n), \m) ® \n — 1)) and neglecting terms of order 
fi/A <C 1. Up to this approximation, the states |m)(g|n— 1) 
(resp. \g) ® \n+ 1)) remain decoupled from the rest for Ui 
(resp. U 2 ). One gets: 



U T |e) =at (e jAt » 



cos 



62 v^ 7 



sin(6>iVNTl) 
2VN + I 



\g) 



( 



e -iA(T-t s )/2 a bl " 2 



N 



COS 



^h) (11) 



with 6<i = VL{T - t s )/2 and 6 2 = Ot s . To make the 
target field state |n) invariant under the Kraus map (8), 
we make it invariant under Ut|c). This is achieved with 
the trapping-like condition: 

Z,7T 7T 

{T-t a ) = —== «• fli = -===. (12) 



For simplicity, we take A t s — nominally 2 . The value of 
A(T — t s ) is irrelevant as it drops out of the Kraus map 
for the field evolution. The value of 82 remains to be fixed. 



2 To this end, the value of A can be slightly tuned by adjusting the 
trap that governs field mode frequency uif. Indeed since A 2> Q, 
a small shift in A allows to sensibly tune A t a for all t s that yield 
non-negligible values of 82 = Qt s . 



4. CONVERGENCE ANALYSIS 

Condition (12) makes the subspace "Hg 4rl+3 ' ) invariant by 
M g , M e and M m . We denote ■p( 4n+3 * > the orthogonal 

projection onto Ti.^ n+3 \ Consider a candidate Lyapunov 
function of the form V(p) = trace(/N p), for some function 
/(n) to be determined. Then (8) and (11) lead to 

V($(p)) - V( P ) = 



trace ^ p cos 2 (^) sin 2 (/3 N ) ( /(N - I) - /(N) ) 

+ trace (p sin 2 (a N ) cos 4 (^) (/(N + I) - /(N) ) ) 

with a N = tt\J^$y and /3 N = 6> 2 VN/2. Note that the last 
line of the above equation vanishes for p = |4n + 3)(4n + 3| 
and any function /n, reflecting the decoupling of T-L^ n+3 ^ 
from the remainder of the Hilbert space. To formally 

restrict ourselves to H^ n+3 \ we can take f(n) = /(4n+3) 
for all n > 4n+3, such that /(N-I) -/(N) and/(N+I)- 
/(N) vanish on H^ +i} . Further take jj G (0, 1), /(n) = 0, 
/(n + 1) = f(n - 1) = 1 and set 

/(n-1) = /(n) + 77 sin 2 ^ cos 2 % (/(n) - /(n+1)) 

for < n < n , 
/(n+1) = f(n) + V sin 2 4f (/(n) - /(n-1)) (13) 

for n < 77 < 4n + 3 . 

Then V($(p)) — V(p) = trace(<7N p) with q n = for n = n, 
n > 4n + 3, and 

g„ = sin 2 a n cos 4 ^ (77 sin 2 %■ - 1) (/(n)-/(n+l)) 

for < 77 < n , (14) 

<&, = sin 2 ,3 n cos 2 ^ (77 sin 2 ^ cos 2 % - 1) (/(n)-/(n-l)) 

for n < 77 < 4n + 3 . 
We then have the following convergence result. 

Theorem 1. Consider the dynamics (8) where M g , M e and 
M m are defined by (11) with (12), and its restriction to 

density operators p with support in Hq 4 ™ +3 ' 1 . Assume that 
02 ^ Wv^for all (n,k) G {1, 4n + 3} x N. Then V(p) 
built with (13) is a strict Lyapunov function: for any po 

with support in 'H^ n+3 \ pk converges towards the fixed 



point p 



\n){n\. 



Proof: Thanks to the assumption on 9 2 , (13) implies 
/(n) > f(n - 1) for 77 < 77 < 477 + 3 and f(n) > f(n + 1) 
for < n < 77. Then the same assumption ensures that q n 
is strictly negative for all 77 < 4n+3 except n = n. Writing 
p = pj I (spectral decomposition with pj > 0, 

^2jPj — 1) then yields 

V(*(p)) - V(p) = pM\9** < 

3 

unless T~'b 4n+3 ' > \'4>j} = Sj \n) with Sj G C for all j, implying 
v (m+3) pV (4n+3) = for somc p e R ^ Thus y i s a 



□ 



strict Lyapunov function in 'Hg 4rl+3 ' ) as stated. 

Remark: The above construction of /(n) can in fact be 
extended to 77 < 9f7 + 8. For generic 82, such that sin(2/3„) 



never vanishes, the resulting V is then a non-strict Lya- 
punov function on the invariant Hilbert subspace 'H^ n+8 " > , 
with the Lasalle invariance principle ensuring convergence 
of p to a poo supported on the subspace spanned by 
\m), |9m + 8). Considering even larger Hilbert spaces, one 
shows that p converges to a state with nonzero population 
only on Fock states \(2l + l) 2 (n + 1) — 1) for I integer. 
The experiments working in very low temperature envi- 
ronments ensure that decoherence usually pulls the field 
state towards vacuum (Fock state |0)) and depopulates 
high-number Fock states. 

Tailored interaction (10) with the atomic stream thus 
makes the otherwise isolated field converge to \n) for 
virtually all 8 2 > 0. This confirms the interest of this 
symmetric product operator approach. We next analyze 
how the latter can be optimized to strengthen convergence 
w.r.t. external disturbances. We therefore add to the 
evolution model a typical relaxation disturbance, also 
called decoherence: evolution between consecutive atomic 
samples becomes 

Pk+1 = $(p fe ) - ^[N$(p fc ) + $(p fe )N - 2a<5(p fe )at] (15) 

-£[(N + I)*(p fe ) + *(p fe )(N + I) - 2at<f>(p fe )a] 

withT + = nn th T s < r~ =n{l+n th )T s < 1. The terms in 
r~ and r + model interaction with a thermal environment 
that induces photon annihilation and creation respectively. 
Parameter k models coupling strength to this environ- 
ment, while nth is the average number of thermal pho- 
tons per mode in the environment. The invariant-subspace 
structure no longer rigorously holds, but we can invoke a 
physical argument to still truncate our computations at 
77 < 9f7 + 8, ensuring that the discretized model remains 
valid and all parameters in the following are well-behaved. 

Writing the operator as a matrix with components 
[Pk]a,b, for a > and b > 0, the evolution (15) couples 
component [pfc] a ,b only to components [pk]a+i,b+i i- e - on 
the same diagonal of the matrix. To analyze how the 
fidelity (n|pfe|n) to target \n)(n\ evolves, we may thus 
reduce our investigation to the vector containing the 
principal diagonal of pk, that is [rk]a = [pk]a.a for a = 
0, 1, 2, ... . Evolution (15) yields 

r k +\ = B ■ A-r k 

with A and B tridiagonal real positive non-symmetric 
matrices, representing reservoir and decoherence respec- 
tively. The upper, lower and principal diagonals of B 
(respectively A) have elements b n = T~ 77 (rcsp. d n = 
sin 2 j3 n cos 2 ^f- ) for 77 G {1,2,...}; c„ = T+ 77 (rcsp. e n = 

sin 2 a n cos 4 4p ) for 77 G {0,1,2,...}; and 1 — c n — b n 
(resp. 1 — d n — e„). A and B are thus column-stochastic, 
reflecting conservation of trace(p) = l T r where 1 is the 
vector of all ones. Therefore B ■ A has at least one eigen- 
value 1. The corresponding 1-cigenvector r* characterizes 
photon populations for a field pointer state under reservoir 
and decoherence. 

We try to approximate it by viewing B ■ A as a small 
perturbation of A. Denote x° the r-vector corresponding 
to \n)(n\, let R° = A-Iand R 1 = (B-T)(A-I) + (B-I). 
Then we know that R°x a = with l T x° = 1, i.e. x° is the 
steady-state solution r* when B = I (no perturbation) . In 



presence of a thermal environment, we get r* = x° + x 1 
solution of 

{R° + R l ){x° + x 1 ) = 0, l T x l =0. (16) 



Proposition 2. Approximating R 1 x 1 rs 0, problem (16) 
can be explicitly solved. It gives 1+a;^ as estimated fidelity 
of poo to our goal, where 



(n\p k \n) 



in-l+l 



, n m 

m=l 1=2 



e-n+l-l 
dn+1 r^O t=2 dfl+l 



En 



An approximation error of order {x\) 2 is expected. 

Proof: Solving (16) approximated (note that dn = &n = 0) 
fixes the components of x , except ir^ which remains free: 



— On / &n~l 
x n + l ~ c n/df L +l 



x n = x n+1 \/n<n~l 

1 e„_i 1 
x„ = —j— x n _ x Vn > n + 1 



□ 



Then (16) determines ir^ < 0. 

Proposition 2 can be used to optimize 02 for maximal 
fidelity in presence of dccohcrcncc, see next Section. Note 
that 8 2 — (the case of Rempe et al. [1990]) would yield 
d n = for all n. The small- 2; 1 approximation would then 
lead to an invalid recurrence, suggesting that decoherence 
leads to large x\ i.e. the method of Rempe et al. [1990] is 
poorly robust. 

5. SIMULATIONS 

For the simulations we consider realistic parameters cor- 
responding to the cavity quantum electrodynamics setup 
at Ecole Normale Superieure, Paris. See e.g. Haroche and 
Raimond [2006], Deleglise and et al. [2008] for detailed 
explanations. Vacuum Rabi frequency is given by f2/27r = 
50 kHz and we take A w 100f2. We numerically compute 
propagators for the exact interaction (7) and u given by 
(10). For each simulation we slightly adapt A to opti- 
mize fidelity, see footnote 2. Evolution of pk is computed 
starting from a vacuum initial state po = |0)(0|. Figure 1 
illustrates the good working principle of the method, de- 
spite our approximate reasoning: pk+i = &(pk) converges 
essentially exactly to p^ = \n)(n\ for all n £ {1,2,..., 8}. 
The simulations are run with an arbitrary 82 = \j\Jn~. 



We next add decoherence. We take 1/k — 0.1 s and cryo- 
genic nth — 0.05, corresponding to current high-standard 
experiments. This environment-induced decoherence is in 
competition with our reservoir strength, that is mainly 
the time between consecutive atoms. The latter can re- 
alistically 3 be set to T s — 60 ps, but there is only a 
0.3 probability that an interacting atom is indeed present 
when trying to send one. The expected evolution from pk 
to Pk+i is thus similar to (a less approximated version 
of) (15) except that operator 0.71 + 0.3 $ replaces 
Figure 2 represents the evolution of the diagonal elements 
of p with the reservoir of Rempe et al. [1990] (top) and 
with our symmetric-interaction reservoir (bottom), for a 

3 Although the required values of Q\ and 82 allow smaller T, we are 
experimentally limited to periods of the order of T B . 
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Fig. 1. Fidelity (n|/9fc|n) of pk to target state \n) as a 
function of atomic interaction fc, for fi — 1,2,..., 8 
(going from top to bottom curve) and r~ = T + = 
i.e. no environmental disturbance. 

target n — 3. For the latter, after the |3) state has built 
up, it progressively drives away towards |15) = \4n + 3) 
and further to higher photon numbers (here accumulating 
where our Hilbert space is truncated). In contrast our 
scheme clearly stabilizes a state close to |3). Achievable 
steady-state fidelities are represented on Fig. 3. We also 
show the values estimated from Proposition 2. They agree 
with simulations within the expected approximation error, 
confirming our theoretical analysis. The selected optimal 
values of 82 seem to satisfy 6*2"\/H s» 37r/4. 

Finally, we consider errors in model parameters as another 
important robustness criterion. Small variations in 82 
appear to give no detectable effect, indicating a rather flat 
optimum in 82- Absolute errors of ±7r/8 on At s also barely 
affect fidelity (less than 1%). The fidelities obtained with 
our reservoir subject to errors in 9\ (and still in presence 
of decoherence) are represented for each n by the 3rd and 
4th bars from left on Fig. 3. Relative errors of ±2% have 
a detectable but tolerable effect; especially, setting 9\ in 
an interval centered slightly below the ideal value seems 
to be a good robustness compromise. The same ±2% error 
on 0\ would be drastically detrimental with the scheme 
of Rempe et al. [1990]: even in absence of a disturbing 
environment (r~ = T + = 0) the fidelity to |n) would 
then quickly decrease towards zero, dropping to about 0.15 
already after 0.1 s and below 0.01 after 0.25 s. Overall our 
construction of atom-field interaction as a product of non- 
commuting transformations thus leads to a significantly 
more robust stabilization scheme. 

6. CONCLUSION 

This paper proposes an 'engineered reservoir' to stabilize 
Fock states of a quantum harmonic oscillator thanks to its 
interaction with a stream of three-level systems. We use a 
single control signal to tailor a time- varying Hamiltonian 
interaction. We prove that the resulting propagator yields 
a Kraus map that robustly 'traps' the spring state on a 
desired Fock state. This seems to be a potentially practical 
open-loop alternative to measurement-based feedback for 
achieving high-fidelity stabilization of Fock states. The 
proposed method admits several variations for experimen- 
tal implementation. Among others, a similar stabilizing 
effect is obtained by using a stream of two-level atoms each 
undergoing one of two different interactions. Given the 
practical possibilities shown in the present and previous 
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Fig. 2. Evolution of the diagonal elements of pk in presence 
of a disturbing environment like in (15), with the 
reservoir of Rcmpe et al. [1990] (top) and with our 
symmetric- interaction reservoir (bottom), for n = 3. 
Graduation is given as a function of time = k ■ 60 ps. 
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3. Fidelity (n\poo\n) of the steady state p^ in presence 
of a disturbing environment to target states \n), for 
n = 1,2,..., 8. The left- most bars represent fidelity 
with our simulated reservoir, the second are the values 
predicted by Proposition 2, and the right- most ones 
are obtained after 4 s with the reservoir of Rempe 
et al. [1990]. The third and fourth bars from left 
correspond to the same conditions as the first ones, 
but incorporating a systematic error of —2% and +2% 
respectively on f?i . 



papers, future work could further investigate systematic 
methods for designing products of interactions that ro- 
bustly stabilize quantum states. 
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